e02dcf

Purpose

Synopsis

[lamda,mu,c,fp,wrk,iwrk,nx,ny,ifail] = e02dcf(x,y,f,s,lamda,mu<,start,wrk,...
iwrk,nx,ny,lwrk,liwrk,ifail>)

Description

 
 This routine determines a smooth bicubic spline approximation 
 s(x,y) to the set of data points (x ,y ,f   ), for q=1,2,...,m  
                                    q  r  q,r                  x
 and r=1,2,...,m .
                y
 
 The spline is given in the B-spline representation
 
                       n -4 n -4                             
                        x    y                               
                       --   --                               
               s(x,y)= >    >   c  M (x)N (y),                 (1)
                       --   --   ij i    j                   
                       i=1  j=1                              
 
 where M (x) and N (y) denote normalised cubic B-splines, the 
        i         j                                          
 former defined on the knots (lambda)  to (lambda)    and the 
                                     i            i+4        
 latter on the knots (mu)  to (mu)   .  
                         j        j+4                          
 
 The total numbers n  and n  of these knots and their values 
                    x      y                                
 (lambda) ,...,(lambda)   and (mu) ,...,(mu)   are chosen 
         1             n          1         n            
                        x                    y           
 automatically by the routine. The knots (lambda) ,...,
                                                 5    
 (lambda)     and (mu) ,...,(mu)     are the interior knots; they 
         n -4         5         n -4                             
          x                      y                               
 divide the approximation domain [x ,x  ]*[y ,y  ] into (
                                   1  m     1  m        
                                      m        m        
 n -7)*(n -7) subpanels [(lambda) ,(lambda)   ]*[(mu) ,(mu)   ], 
  x      y                       i         i+1       j     j+1  
 for i=4,5,...,n -4, j=4,5,...,n -4. Then, much as in the curve 
                x               y                              
 case (see E02BEF), the coefficients c   are determined as the 
                                      ij                      
 solution of the following constrained minimization problem:
 
 minimize
 
                            (eta),                             (2)
 
 subject to the constraint
 
                        m   m                                
                         x   y                               
                        --  --          2                    
               (theta)= >   >  (epsilon)   <=S,                (3)
                        --  --          q,r                  
                        q=1 r=1                              
 
 where  (eta) is a measure of the (lack of) smoothness of s(x,y). 
              Its value depends on the discontinuity jumps in 
              s(x,y) across the boundaries of the subpanels. It is
              zero only when there are no discontinuities and is 
              positive otherwise, increasing with the size of the 
              jumps.
 
        (epsilon)    denotes the residual f   -s(x ,y ),
                 q,r                       q,r    q  r 
 
 and    S    is a non-negative number to be specified by the user.
 
 By means of the parameter S, 'the smoothing factor', the user 
 will then control the balance between smoothness and closeness of
 fit, as measured by the sum of squares of residuals in (3). If S 
 is too large, the spline will be too smooth and signal will be 
 lost (underfit); if S is too small, the spline will pick up too 
 much noise (overfit). In the extreme cases the routine will 
 return an interpolating spline ((theta)=0) if S is set to zero, 
 and the least-squares bicubic polynomial ((eta)=0) if S is set 
 very large. Experimenting with S-values between these two 
 extremes should result in a good compromise. 
 
 Values of the computed spline can subsequently be computed by 
 calling E02DEF or E02DFF.
 

Parameters

e02dcf

Required Input Arguments:

x (:)                                 real
y (:)                                 real
f (:)                                 real
s                                     real
lamda (:)                             real
mu (:)                                real

Optional Input Arguments:                       <Default>

start (1)                             string   'c'
wrk (lwrk)                            real     if lower(start)=='c';wrk=...
                                        ...    zeros(lwrk,1);end
iwrk (liwrk)                          integer  if lower(start)=='c';iwrk=...
                                        ...    zeros(liwrk,1);end
nx                                    integer  if lower(start)=='c';nx=0;end
ny                                    integer  if lower(start)=='c';ny=0;end
lwrk                                  integer  e02dcf17(length(x),...
                                        ...    length(y),length(lamda),...
                                        ...    length(mu))
liwrk                                 integer  e02dcf19(length(x),...
                                        ...    length(y),length(lamda),...
                                        ...    length(mu))
ifail                                 integer  -1

Output Arguments:

lamda (:)                             real
mu (:)                                real
c (:)                                 real
fp                                    real
wrk (lwrk)                            real
iwrk (liwrk)                          integer
nx                                    integer
ny                                    integer
ifail                                 integer